The Interacting Branching Process as a Simple Model of Innovation 
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We describe innovation in terms of a generalized branching process. Each new invention pairs 
with any existing one to produce a number of offspring, which is Poisson distributed with mean p. 
Existing inventions die with probability p/r at each generation. In contrast to mean field results, 
no phase transition occurs; the chance for survival is finite for all p > 0. For t = oo, surviving 

processes exhibit a bottleneck before exploding super-exponentially - a growth consistent with a 
law of accelerating returns. This behavior persists for finite r. We analyze, in detail, the asymptotic 
behavior as p ^ 0. 
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Networks of elements creating new elements can be 
found in diverse conditions - from the origin of life to 
artistic expression. As such, innovation and discovery are 
intrinsic to life as well as human experience. Economic, 
technological or social innovations include introducing a 
new good, method of production, form of governance, 
etc. Clearly, innovations do not happen in a vacuum. 
They form contingent and interconnected webs where ex- 
isting discoveries foster the creation of new ones - leading 
to a network of self-perpetuating, autocatalytic activity. 
Landmark innovations ignite radical change - avalanches 
of new discoveries that were previously unthinkable [1]. 
These include, for instance, the World Wide Web, or the 
bursts of creativity associated with the Renaissance. But 
innovation is not restricted to human society. It is an 
archetype for any co-evolutionary dynamics. Indeed the 
intermittent, bursty pattern of innovation in human his- 
tory resembles punctuated equilibrium observed in bio- 
logical evolution - where episodes, such as the Cambrian 
explosion, dominate the history of life's diversity [2]. 

Previous theoretical approaches to this phenomenon 
have concentrated on mean field models, where explo- 
sions of innovation were attributed to a phase transi- 
tion [3-5]. This transition separates a regime where ac- 
tivity always dies out after giving rise to only a few el- 
ements from a phase where large numbers of elements 
(or inventions) can come forth. The phase transition 
depends both on the number of initial - or primeval - 
inventions and the "mating" probability to create new 
elements from interactions between existing ones. 

Here we take a different approach that explicitly con- 
siders fluctuations in a microscopic model of innovation. 
Pairs of inventions mate to create new ones. The relevant 
criterion is the chance for the process to survive forever, 
or, in other words to escape extinction. In contrast to 
conclusions based on mean field arguments, our analy- 
sis establishes that fluctuations always lead to a finite 
probability to escape extinction for any mating probabil- 
ity p ^ 0, and any positive number of initial inventions. 
Hence no phase transition occurs. For small p, the pop- 
ulation exhibits a long bottleneck, where it grows slowly 



with time. But after this quiescent period, surviving pop- 
ulations explode faster than exponentially. Indeed both 
human populations [6] and technological advances (such 
as computing speed for fixed cost [7]) are known to ex- 
hibit a 'law of accelerating return' [8] - where growth 
occurs faster than exponentially. We associate the sud- 
den proliferation of innovations following a landmark in- 
vention with such a super-explosion, rather than with a 
phase transition. We derive analytic results, including 
exact scaling laws in the limit p — > 0, for the bottle- 
neck and transition to super-explosive growth that are 
confirmed by numerical simulations. 

At its most basic level, innovation can be represented 
as a growing phylogenetic network. Phylogenetic net- 
works [9] are generalizations of trees to the case of more 
than one ancestor. Such trees find wide uses to depict 
evolutionary relationships between e.g. genes or species, 
cultures [10], or languages [12]. Each node in the network 
represents an element in the population of innovations. It 
has one or more parent nodes, and can branch to produce 
daughter nodes (new innovations) at future generations. 

Independent branching of each node corresponds to the 
Galton- Watson (GW) branching process. Assume that 
starting with a single root node, each node on the tree 
independently produces a number of daughter nodes that 
is Poisson distributed with mean fi. The whole tree goes 
extinct only if each of the subtrees that follow the root's 
daughters dies. Hence, the survival probability Z satis- 
fies 1 — Z = e~^^ , with a non-zero solution only when 
H > 1. This phase transition that the GW process ex- 
hibits at /U = 1 is a general property of growth processes 
with independent branching. 

In contrast, we recognize that innovation is a histor- 
ically contingent process driven by interactions. Hence 
we put forward a more relevant process - the interacting 
branching process (IBP), as a prototypical model. The 
total population at generation g is Ng. It consists of sur- 
vivors from previous generations and new nodes created 
in generation g. At update -|- 1, if Sg+i "young" nodes 
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are born, and Kg^i "old" nodes are killed, 

Ng + ^=Ng+Sg+^-Kg+,. (l) 

Birth happens because each of the Sg nodes in genera- 
tion g can mate with each of the Ng nodes currently alive 
(including young and old ones, and even itself). Each 
mating pair produces a Poisson distributed number of 
offspring with mean p. This makes the total number of 
offspring of a node in a generation g Poisson distributed 
with mean Ug = pNg. Hence the size of the next genera- 
tion Sg+i (given Sg and Ng) is also Poisson distributed, 
with mean 5*^+1 = pNgSg. Note that Ug in the IBP 
depends on Ng, unlike the GW process. 

Killing nodes incorporates the possibility that inven- 
tions become obsolete. Killing also happens stochasti- 
cally: After creating the new (g + l)*'* generation, we kill 
and remove each of the Ng old individuals, with proba- 
bility p/t. Hence the average Kg+i = {p/T)Ng. 

In a mean field analysis, the stochastic variables Sg and 
Kg are replaced by their average values Sg and Kg . If r = 
p no nodes from older generations survive. Then Ng = Sg 
for all g and -Sg+i = pS"^. If the initial value < 
the process always dies out, while it escapes extinction if 
> 1/p - indicating a phase transition. However, due to 
fluctuations, some populations can survive for any finite 
p, even when Sq = \. The probability to survive goes to 
zero when p — > 0, but it is non-zero for any p > Q. In 
addition, there exists a critical value of r such that for 
T > Tc the survival probability of the IBP has a precisely 
known essential singularity as p 0. We also find that 
fluctuations drive super-explosions for t < Tc] while for 
T > Tc fluctuations are irrelevant. The precise value of 
Tc depends on (tc ~ 0.391 for = 1), but it is finite 
as long as So is finite. For brevity, we restrict So = 1 in 
what follows. 

For T = oo innovations never die. Kg = for all g, 
and Ng — Yl^g'={) ^g'- Before the occurence of an empty 
generation the branching ratio Ug > pg is an increas- 
ing function of g with average slope > p. Hence, for 
some long lived processes Ug must exceed one at finite 
g. Comparison with the GW process shows that after 
Ug has irreversibly exceeded unity, the population has a 
non-zero probability to escape extinction. Notice that 
this is a rigorous result and contradicts the conchisions 
of [4, 5]. Unlike the GW process, the average IBP popu- 
lation size does not explode exponentially with time. It 
super-explodes: Soon after Ug exceeds one, the popula- 
tion size grows faster than exponentially. In fact this 
is the almost certain fate of an IBP that starts with 
So = No ^ 1/p. Then we can ignore fiuctuations and 

set Sg+l =pNgSg, OT 




indicating that pSg and hence Sg grows faster than 
(pSo)^^ - or faster than exponentially. We say that the 



IBP reaches the super-explosive phase when Sg increases 
faster than exponentially with g, and now study how the 
IBP approaches this regime. 

For any t, if a population has survived to generation g, 
it dies at the {g + 1)*'' generation if none of the Sg nodes 
have any off'spring, an event that happens with probabil- 
ity e~P^9'^3. Thus if Z{g) is the probability to survive 
to the g*'' generation, Z{g + 1) = [l - (g-P^s'^s)') z{g). 
The prime indicates that the average is restricted to pop- 
ulations that survive to g. These are "conditioned-to- 
grow" populations, which do not include empty genera- 
tions. Iterating this expression gives 

Zig)=liil-{e-^''^''^')') . (2) 

Setting (? = 00, we get the probability that the process 
escapes extinction, Z = Z{oo). 

For p 1, almost all populations die out after a few 
generations, frustrating a direct numerical approach to 
obtain Z. However, Eq. (2) motivates a computationally 
efiicient method that conditions on surviving populations 
and ignores those that die out. The conditioned-to-grow 
populations entering the expectation value in Eq. (2) 
have a Poisson offspring distribution 

p-pNgSg I Q \m 

Prob'[5,+, = m\NgSg] = ^ 'J^ , (3) 

for m > 1. The trunctation at m = 1 corresponds to 
conditioning on survival. We record Ng, Sg for each gen- 
eration to obtain the death probability {e~^^^^^)'. We 
now define a rescaled time t = pg, and write the death 
probability as (^e^P^s'^gy = e^*2?(t), in which case the 
logarithm of the survival probability Z{t/p) obeys 

t/p 

limp log Z(t/p) = limp y ln(l - e-^^s'^f')' 

S'=0 

= ^ dt'ln (l - e-*'p(i')) , (4) 

where the last expression is manifestly independent of p. 

Fig. 1 shows the scaled logarithm of the survival prob- 
ability as a function of rescaled time for r = oo. It ex- 
hibits near-perfect collapse for different p -C 1, in agree- 
ment with Eq. (4). Fig. 1 also indicates that the scaled 
death probability goes to zero at t = ~ 0.67 for p 0, 
and concomitantly the survival probability becomes con- 
stant for all subsequent times. Thus ^ ~ e~^/^ for p — )• 0, 
with c « 0.89. These numerical results support our claim 
that for r = 00 the IBP has a finite probability to escape 
extinction for any p ^ 0. 

We now show how the branching ratio Ug = pNg 
approaches unity from below, before the IBP super- 
explodes. Multiplying Eq. (3) by Proh'[Sg\Ng] and sum- 
ming over Sg gives Prob' [5*^+1 = m\Ng] conditioned on 
growth for m > 1 and fixed Ng. Using continuous time t 
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FIG. 1. Death and survival probabilities, (color online) 
Diamonds: rescaled rate T>{t) at which surviving populations 
die vs. scaled time t = pg, for two different values of p ^ 1, 
W~^{blue), and 10~^(red). As p — >■ 0, the death probability 
converges to a limiting function that vanishes a.t t — te ~ 0.67. 
Circles: —p log Z{t/p) vs. t, for the same two values of p. 



and changing notation, Eq. (3) leads to 

,^U{t)['3 ,1-)- 

S' 



Puit){S,t + dt) = Y,Pu(t){S\t)- 



p-S'U(t) 



Setting the conditioning variable U (t) = U, numerical it- 
erations of this equation quickly converge to a stationary 
distribution Pu{S) for any U < 1. Since U{t) involves an 
integral of S(t), it has both lower fluctuations and slower 
variations compared to S, as long as U{t) < 1. Hence, 
by the law of large numbers, U (t) can be replaced by its 
mean U'{t) over different realizations surviving to time t. 
Then the distribution of S for a given u = U' is obtained 
from the stationary solution of 



-S' u 



[s'uY 



o — S'u 



(5) 



Noting that as neither p, time, nor the age of the nodes 
enter Eq. (5), we expect it to be valid for all p ^ 0, and 
also for all r - as long as u < 1. 

We have tested the distribution given by Eq. (5) 
against averages over 10^ realizations of conditioned-to- 
grow IBP populations for small values of p. Fig. 2 
shows our result for the mean s = S' = J2s' ^'Pu{S') vs. 
u. The agreement is excellent not only for the example 
T = oo shown, but also for all other t > Tc ^ 0.391. The 
latter condition is explained next. 

For finite r, the branching ratio evolves a,s U — S{t) — 
K{t), where K{t) is the number of nodes killed at time 
t. Neglecting fluctuations gives ii — v{u) where v{u) = 
s{u) — u/t. Fig. 2 shows that v{u) > for all u, if 
T > Tc « 0.391. As a result, u{t) increases and eventually 
super-explodes, so populations have a finite probability 
to escape extinction. 

For r < Tc, u/t intersects s{u) at two values of u: ui 
and U2. The smaller of these, ui, gives the maximal value 
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FIG. 2. Mean generation size as a function of branch- 
ing ratio, (color online) The solid black curve is the numer- 
ical prediction for the mean generation size s vs. the mean 
branching ratio u (both conditioned on growth) as obtained 
from Eq. (5). The diamonds are results of numerical simula- 
tions of the IBP. We used p — 10~® and r — oo, and averaged 
over 10^ conditioned-to-grow populations. Agreement is ex- 
cellent. Also plotted are three lines showing u/t on the y-axis, 
for three values of r from top to bottom: (1) 0.16, sub-critical 
T < Tc, (2) critical t — Tc = 0.391, and (3) 0.5, super-critical 
T > Tc, respectively. 




FIG. 3. Mean time te to enter the super-explosive 
phase, (color online) The red circles (blue diamonds ) are 
the results of numerical simulations at p = 10~® (10~^) for 
the mean time to enter the super-explosive phase, the black 
curve going through them a fit oc (t — r^)"^'''^. 



of u{t) reached before the process goes extinct. Including 
fluctuations in a Langevin approach, ii = v{u)+p2U!(u)(^, 
with C a standard Gaussian noise [13]. Now the state 
with u = ui is met ast able: For any flnite p, a flnite 
fluctuation can kick the system over the potential barrier 
to the unstable value U2. Beyond U2, s is larger than u/t 
for all M, and surviving processes again super-explode. 

This scenario is supported by results of numerical sim- 
ulations, which we present in Fig. 3. The time, te, to 
super-explode diverges as (r — Tc)^^/^. It is indepen- 
dent of p, for T > Tc provided (r — Tc) is sufficiently 
large. This divergence is obtained analytically by ex- 
panding u about Umin ~ the value where v{u) reaches 
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FIG. 4. Mean first passage (generation) times (MFPT) 
to reacli a generation of size m, for r = oo. MFPTs 
for m = 2, 3, 4, 5 are shown from bottom to top respectively, 
along with the predicted values from Eq. (7). For the range 
of p studied, numerical results obtained by averaging over 10® 
surviving processes agree with predictions for m = 2,3. For 
m = 4, 5 they converge to the predicted values as p — )■ 0. 



its minimum Vmin- Writing u = Umin + e, one gets 
de/di = 1 + ^2/2 + o(?2). Here, i = e^^^Aw, 

i = ty^s'^i^Vrntn, and s^j„ = d'^s/du'^\u=u^i„- Surviv- 
ing populations spend most of the time around Umin be- 
fore reaching u = 1 and super-exploding. Hence, the 
time to explode is te ^ (Smm^m™)""^^^- At t = Tc, 
Vmin = 0. A Taylor expansion near gives ^ t—Tc, 
and hence ^ {t — t^~^I'^ , in agreement with Fig. 3. 
For finite p, the scaling breaks down when the drift 
term Vmin becomes comparable to the noise, or when 



For T < Tc, the mean time to enter the explosive phase, 
te, is no longer independent of p in the limit p — 0. For 
e.g. T = p, one has Ug = pSg, and U > 1 requires 
Sg > 1/p. The next paragraph shows that the average 
time to first reach Sg > 1/p is te ~ r(l/p)p2-(i/p). 

Finally we estimate the mean first passage (rescaled) 
time, Tm{p) = pGmip), to a generation of size m or larger 
in populations conditioned to reach such a generation 
size. We derive upper bounds by considering only the 
single most probable path of evolution, which become 
exact in the limit p — ?> 0. For any r and small p, the 
most likely surviving process before super-exploding is a 
simple chain where Sg = 1 for all g. For populations 
conditioned to reach Sg > m for some 5, the most likely 
shape is a chain up to 5 — 1, and a fan-out from Sg-i =^ 1 
to S„ — m during the last generation. All other shapes 



are reduced by factors of p. Hence for t = p (so Ugi = p 
for all g' < g in the chain) this gives the same relative 
chance p"^^^/m\ to reach Sg > m { compared to Sg = I 
) at any g. The probability ^™ that Sg > m has not been 
reached evolves as = (^1 — p'"~^/m!) and thus 



Qm{p) 



Tmjp) 
P 



m! p^ ™ (for T = p) 



(6) 



The time can be obtained by setting m = I /p. 

For T = 00, Ug ~ pg is no longer independent of g. The 
relative chance to reach Sg > m { compared to Sg ^ 1 
) is (TO)'"-Vm!, and Cgli « (1 - (P5)™" 7^!)^™. Inte- 
grating over g reveals the (7-dependence of and that 



Gmip) oc 



1-1/m 



(for T = 00) 



(7) 



This is compared in Fig. 4 to direct simulations of 
the IBP. It describes the numerical results perfectly for 
p — 0, and gives upper bounds for finite p, as it should. 
Eq. (7) clarifies how surviving populations evolve in the 
limit p — > 0. First, only one individual is born in each 
generation. Generations of size two start appearing after 
p-1/2 generations, followed by the first appearance of a 
generation of size three after p~'^/^ generations and so on. 
Finally after 1/p generations, generations of size s — )• 00 
appear, indicating the onset of super-explosive growth. 

We have described autocatalytic networks of innova- 
tion in terms of an interacting branching process (IBP). 
In contrast to standard branching processes, two parents 
are needed to generate offspring. The IBP is both an- 
alytically and numerically tractable. In mean field the- 
ory, it shows a phase transition, which disappears due to 
fluctuations in a rigorous microscopic treatment. When 
the probability p to make new innovations from any two 
existing ones is vanishingly small, we find universal be- 
havior that is independent of p. A super-explosive phase, 
where the rate of new inventions grows faster than expo- 
nentially follows a long quiescent bottleneck for p <C 1. 
This dynamics resembles the Dark Ages preceding the 
Renaissance or the quiescent times between bursts of spe- 
ciation in punctuated biological evolution. We speculate 
that our models unfolds large scale properties of any co- 
evolutionary dynamics. Indeed its super-explosive be- 
havior is consistent with the law of accelerating returns 
found in technological progress [7, 8]. 
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